GEOG_4

Trajectory Analysis

1. Data Preparation

# Read data, I choose African Lions as the animal
animal <- read.csv('~/Downloads/African Lion - Panthera Lion Central Kalahari Game Reserve.csv') # name the file as animal

head(animal)

##      event.id visible               timestamp location.long location.lat
## 1 33042889052    true 2008-12-15 05:00:00.000      23.65028    -21.22511
## 2 33042889053    true 2008-12-15 06:00:00.000      23.65590    -21.22162
## 3 33042889054    true 2008-12-15 07:00:00.000      23.65950    -21.22039
## 4 33042889055    true 2008-12-15 08:00:00.000      23.65954    -21.22042
## 5 33042889056    true 2008-12-15 09:00:00.000      23.65981    -21.22146
## 6 33042889057    true 2008-12-15 10:00:00.000      23.65984    -21.22146
##   eobs.temperature gps.dop height.raw sensor.type
## 1               29     1.8     987.00         gps
## 2               30     2.0     983.90         gps
## 3               31     1.8     978.80         gps
## 4               33     3.4     990.30         gps
## 5               35     2.2     987.63         gps
## 6               36     3.6     987.52         gps
##   individual.taxon.canonical.name tag.local.identifier
## 1                    Panthera leo             GSM05751
## 2                    Panthera leo             GSM05751
## 3                    Panthera leo             GSM05751
## 4                    Panthera leo             GSM05751
## 5                    Panthera leo             GSM05751
## 6                    Panthera leo             GSM05751
##   individual.local.identifier
## 1                        1003
## 2                        1003
## 3                        1003
## 4                        1003
## 5                        1003
## 6                        1003
##                                                   study.name
## 1 African Lion - Panthera Lion Central Kalahari Game Reserve
## 2 African Lion - Panthera Lion Central Kalahari Game Reserve
## 3 African Lion - Panthera Lion Central Kalahari Game Reserve
## 4 African Lion - Panthera Lion Central Kalahari Game Reserve
## 5 African Lion - Panthera Lion Central Kalahari Game Reserve
## 6 African Lion - Panthera Lion Central Kalahari Game Reserve

# check what variable are available
names(animal)

##  [1] "event.id"                        "visible"                        
##  [3] "timestamp"                       "location.long"                  
##  [5] "location.lat"                    "eobs.temperature"               
##  [7] "gps.dop"                         "height.raw"                     
##  [9] "sensor.type"                     "individual.taxon.canonical.name"
## [11] "tag.local.identifier"            "individual.local.identifier"    
## [13] "study.name"

# Unique individuals
individuals <- unique(animal$individual.local.identifier)
individuals

##  [1] 1003 1001 1002 1013 1008 1012 1006 1007 1005 1004 1010 1009 1011

We have 13 individuals in the data.

2. Trajectory Map (static maps)

animal_sf <- st_as_sf(animal, coords = c("location.long", "location.lat"), crs = 4326)

# Switch to interactive mode
tmap_mode("plot") # for your final maps you can use static maps

## ℹ tmap mode set to "plot".

# Basic map with OpenStreetMap + African Lions points
tm_shape(animal_sf) +
  tm_basemap(server = "OpenStreetMap") +
  tm_symbols(
    col = "individual.local.identifier",
    palette = "Set1",
    size = 0.5,
    border.col = "black",     # set outline colour
    border.lwd = 1            # set outline thickness (1 is default)
  ) +
  tm_view(set.view = c(-22, 17, 10))  # Center + zoom

## 
## ── tmap v3 code detected ───────────────────────────────────────────────────────
## [v3->v4] `tm_symbols()`: migrate the argument(s) related to the scale of the
## visual variable `fill` namely 'palette' (rename to 'values') to fill.scale =
## tm_scale(<HERE>).[v3->v4] `symbols()`: use 'fill' for the fill color of polygons/symbols
## (instead of 'col'), and 'col' for the outlines (instead of 'border.col').[v3->v4] `tm_view()`: use set_view instead of set.view[cols4all] color palettes: use palettes from the R package cols4all. Run
## `cols4all::c4a_gui()` to explore them. The old palette name "Set1" is named
## "brewer.set1"

3. Segment Duration, Length and Speed

Prepare the data:

# Check which rows have NA values in coordinates
NAPoints <- which(is.na(animal$location.lat) | is.na(animal$location.long))

# First we need to convert timestamp from factor to datetime
animal$tstamp <- strptime(as.character(animal$timestamp), "%Y-%m-%d %H:%M:%S")

# Then we separate timestamp into individual components
animal$date <- as.Date(animal$tstamp)
animal$year <- year(animal$tstamp)
animal$month <- month(animal$tstamp)
animal$day <- day(animal$tstamp)
animal$hour <- hour(animal$tstamp)
animal$min <- minute(animal$tstamp)
animal$sec <- second(animal$tstamp)

head(animal)

##      event.id visible               timestamp location.long location.lat
## 1 33042889052    true 2008-12-15 05:00:00.000      23.65028    -21.22511
## 2 33042889053    true 2008-12-15 06:00:00.000      23.65590    -21.22162
## 3 33042889054    true 2008-12-15 07:00:00.000      23.65950    -21.22039
## 4 33042889055    true 2008-12-15 08:00:00.000      23.65954    -21.22042
## 5 33042889056    true 2008-12-15 09:00:00.000      23.65981    -21.22146
## 6 33042889057    true 2008-12-15 10:00:00.000      23.65984    -21.22146
##   eobs.temperature gps.dop height.raw sensor.type
## 1               29     1.8     987.00         gps
## 2               30     2.0     983.90         gps
## 3               31     1.8     978.80         gps
## 4               33     3.4     990.30         gps
## 5               35     2.2     987.63         gps
## 6               36     3.6     987.52         gps
##   individual.taxon.canonical.name tag.local.identifier
## 1                    Panthera leo             GSM05751
## 2                    Panthera leo             GSM05751
## 3                    Panthera leo             GSM05751
## 4                    Panthera leo             GSM05751
## 5                    Panthera leo             GSM05751
## 6                    Panthera leo             GSM05751
##   individual.local.identifier
## 1                        1003
## 2                        1003
## 3                        1003
## 4                        1003
## 5                        1003
## 6                        1003
##                                                   study.name
## 1 African Lion - Panthera Lion Central Kalahari Game Reserve
## 2 African Lion - Panthera Lion Central Kalahari Game Reserve
## 3 African Lion - Panthera Lion Central Kalahari Game Reserve
## 4 African Lion - Panthera Lion Central Kalahari Game Reserve
## 5 African Lion - Panthera Lion Central Kalahari Game Reserve
## 6 African Lion - Panthera Lion Central Kalahari Game Reserve
##                tstamp       date year month day hour min sec
## 1 2008-12-15 05:00:00 2008-12-15 2008    12  15    5   0   0
## 2 2008-12-15 06:00:00 2008-12-15 2008    12  15    6   0   0
## 3 2008-12-15 07:00:00 2008-12-15 2008    12  15    7   0   0
## 4 2008-12-15 08:00:00 2008-12-15 2008    12  15    8   0   0
## 5 2008-12-15 09:00:00 2008-12-15 2008    12  15    9   0   0
## 6 2008-12-15 10:00:00 2008-12-15 2008    12  15   10   0   0

# What years do we have?
years <- unique(animal$year)

#Create a birdyear id for so that trajectories of different years can be separated if we so wish.
animal$idyear<-paste0(animal$individual.local.identifier, year(animal$date))

# Get individuals by year
individualYear <- unique(animal$idyear)
individualYear

##  [1] "10032008" "10032009" "10012009" "10012010" "10012011" "10022011"
##  [7] "10022012" "10132010" "10132011" "10132012" "10082010" "10082011"
## [13] "10122010" "10122011" "10062010" "10062011" "10072011" "10052011"
## [19] "10052012" "10042011" "10042012" "10022009" "10022010" "10102009"
## [25] "10122009" "10092009" "10092010" "10112009" "10112010" "10072010"
## [31] "10042010" "10092011" "10112011"

The time period from 2008 to 2011 in the data.

# Convert timestamp to POSIX format
animal$POSIX <- as.POSIXct(animal$timestamp, format="%Y-%m-%d %H:%M:%S")

# Are any times duplicated?
duplics <- which(duplicated(animal$POSIX)==TRUE)

# Yes, we have a few points like this, so let's remove them.
animal <- animal[-duplics,]

# Convert the ducks data frame into an sf object by specifying which columns are coordinates
# Also we need to specify coordinate system, which in our case is WGS1984 (EPSG:4326 - for EPSG codes we just need the number))
animalSF <- st_as_sf(animal, coords=c('location.long','location.lat'), crs=4326)
animalSF

## Simple feature collection with 41614 features and 21 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: 22.37396 ymin: -21.90794 xmax: 23.95741 ymax: -20.84132
## Geodetic CRS:  WGS 84
## First 10 features:
##       event.id visible               timestamp eobs.temperature gps.dop
## 1  33042889052    true 2008-12-15 05:00:00.000               29     1.8
## 2  33042889053    true 2008-12-15 06:00:00.000               30     2.0
## 3  33042889054    true 2008-12-15 07:00:00.000               31     1.8
## 4  33042889055    true 2008-12-15 08:00:00.000               33     3.4
## 5  33042889056    true 2008-12-15 09:00:00.000               35     2.2
## 6  33042889057    true 2008-12-15 10:00:00.000               36     3.6
## 7  33042889058    true 2008-12-15 12:00:00.000               34     2.2
## 8  33042889059    true 2008-12-15 13:00:00.000               34     2.2
## 9  33042889060    true 2008-12-15 14:00:00.000               34     2.0
## 10 33042889061    true 2008-12-15 16:00:00.000               33     1.6
##    height.raw sensor.type individual.taxon.canonical.name tag.local.identifier
## 1      987.00         gps                    Panthera leo             GSM05751
## 2      983.90         gps                    Panthera leo             GSM05751
## 3      978.80         gps                    Panthera leo             GSM05751
## 4      990.30         gps                    Panthera leo             GSM05751
## 5      987.63         gps                    Panthera leo             GSM05751
## 6      987.52         gps                    Panthera leo             GSM05751
## 7      986.41         gps                    Panthera leo             GSM05751
## 8      986.90         gps                    Panthera leo             GSM05751
## 9      986.08         gps                    Panthera leo             GSM05751
## 10     985.46         gps                    Panthera leo             GSM05751
##    individual.local.identifier
## 1                         1003
## 2                         1003
## 3                         1003
## 4                         1003
## 5                         1003
## 6                         1003
## 7                         1003
## 8                         1003
## 9                         1003
## 10                        1003
##                                                    study.name
## 1  African Lion - Panthera Lion Central Kalahari Game Reserve
## 2  African Lion - Panthera Lion Central Kalahari Game Reserve
## 3  African Lion - Panthera Lion Central Kalahari Game Reserve
## 4  African Lion - Panthera Lion Central Kalahari Game Reserve
## 5  African Lion - Panthera Lion Central Kalahari Game Reserve
## 6  African Lion - Panthera Lion Central Kalahari Game Reserve
## 7  African Lion - Panthera Lion Central Kalahari Game Reserve
## 8  African Lion - Panthera Lion Central Kalahari Game Reserve
## 9  African Lion - Panthera Lion Central Kalahari Game Reserve
## 10 African Lion - Panthera Lion Central Kalahari Game Reserve
##                 tstamp       date year month day hour min sec   idyear
## 1  2008-12-15 05:00:00 2008-12-15 2008    12  15    5   0   0 10032008
## 2  2008-12-15 06:00:00 2008-12-15 2008    12  15    6   0   0 10032008
## 3  2008-12-15 07:00:00 2008-12-15 2008    12  15    7   0   0 10032008
## 4  2008-12-15 08:00:00 2008-12-15 2008    12  15    8   0   0 10032008
## 5  2008-12-15 09:00:00 2008-12-15 2008    12  15    9   0   0 10032008
## 6  2008-12-15 10:00:00 2008-12-15 2008    12  15   10   0   0 10032008
## 7  2008-12-15 12:00:00 2008-12-15 2008    12  15   12   0   0 10032008
## 8  2008-12-15 13:00:00 2008-12-15 2008    12  15   13   0   0 10032008
## 9  2008-12-15 14:00:00 2008-12-15 2008    12  15   14   0   0 10032008
## 10 2008-12-15 16:00:00 2008-12-15 2008    12  15   16   0   0 10032008
##                  POSIX                   geometry
## 1  2008-12-15 05:00:00 POINT (23.65028 -21.22511)
## 2  2008-12-15 06:00:00  POINT (23.6559 -21.22162)
## 3  2008-12-15 07:00:00  POINT (23.6595 -21.22039)
## 4  2008-12-15 08:00:00 POINT (23.65954 -21.22042)
## 5  2008-12-15 09:00:00 POINT (23.65981 -21.22146)
## 6  2008-12-15 10:00:00 POINT (23.65984 -21.22146)
## 7  2008-12-15 12:00:00 POINT (23.65982 -21.22146)
## 8  2008-12-15 13:00:00 POINT (23.65981 -21.22145)
## 9  2008-12-15 14:00:00 POINT (23.65987 -21.22147)
## 10 2008-12-15 16:00:00 POINT (23.65985 -21.22145)

# Project data
# Define new projection
crsnew <- st_crs("ESRI:102014")

# Transform data into the new coordinate system
animalSF_proj <- st_transform(animalSF,crs=crsnew)
head(animalSF_proj)

## Simple feature collection with 6 features and 21 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: 2891248 ymin: -7700646 xmax: 2893092 ymax: -7699222
## Projected CRS: Europe_Lambert_Conformal_Conic
##      event.id visible               timestamp eobs.temperature gps.dop
## 1 33042889052    true 2008-12-15 05:00:00.000               29     1.8
## 2 33042889053    true 2008-12-15 06:00:00.000               30     2.0
## 3 33042889054    true 2008-12-15 07:00:00.000               31     1.8
## 4 33042889055    true 2008-12-15 08:00:00.000               33     3.4
## 5 33042889056    true 2008-12-15 09:00:00.000               35     2.2
## 6 33042889057    true 2008-12-15 10:00:00.000               36     3.6
##   height.raw sensor.type individual.taxon.canonical.name tag.local.identifier
## 1     987.00         gps                    Panthera leo             GSM05751
## 2     983.90         gps                    Panthera leo             GSM05751
## 3     978.80         gps                    Panthera leo             GSM05751
## 4     990.30         gps                    Panthera leo             GSM05751
## 5     987.63         gps                    Panthera leo             GSM05751
## 6     987.52         gps                    Panthera leo             GSM05751
##   individual.local.identifier
## 1                        1003
## 2                        1003
## 3                        1003
## 4                        1003
## 5                        1003
## 6                        1003
##                                                   study.name
## 1 African Lion - Panthera Lion Central Kalahari Game Reserve
## 2 African Lion - Panthera Lion Central Kalahari Game Reserve
## 3 African Lion - Panthera Lion Central Kalahari Game Reserve
## 4 African Lion - Panthera Lion Central Kalahari Game Reserve
## 5 African Lion - Panthera Lion Central Kalahari Game Reserve
## 6 African Lion - Panthera Lion Central Kalahari Game Reserve
##                tstamp       date year month day hour min sec   idyear
## 1 2008-12-15 05:00:00 2008-12-15 2008    12  15    5   0   0 10032008
## 2 2008-12-15 06:00:00 2008-12-15 2008    12  15    6   0   0 10032008
## 3 2008-12-15 07:00:00 2008-12-15 2008    12  15    7   0   0 10032008
## 4 2008-12-15 08:00:00 2008-12-15 2008    12  15    8   0   0 10032008
## 5 2008-12-15 09:00:00 2008-12-15 2008    12  15    9   0   0 10032008
## 6 2008-12-15 10:00:00 2008-12-15 2008    12  15   10   0   0 10032008
##                 POSIX                 geometry
## 1 2008-12-15 05:00:00 POINT (2891248 -7700646)
## 2 2008-12-15 06:00:00 POINT (2892276 -7699642)
## 3 2008-12-15 07:00:00 POINT (2892975 -7699222)
## 4 2008-12-15 08:00:00 POINT (2892986 -7699227)
## 5 2008-12-15 09:00:00 POINT (2893087 -7699449)
## 6 2008-12-15 10:00:00 POINT (2893092 -7699448)

# Read X an Y from geometry column as separate coordinates
animalSF_proj$X <- st_coordinates(animalSF_proj)[,1]  
animalSF_proj$Y <- st_coordinates(animalSF_proj)[,2]

# Drop geometry column
animalSF_proj <- st_drop_geometry(animalSF_proj)
head(animalSF_proj)

##      event.id visible               timestamp eobs.temperature gps.dop
## 1 33042889052    true 2008-12-15 05:00:00.000               29     1.8
## 2 33042889053    true 2008-12-15 06:00:00.000               30     2.0
## 3 33042889054    true 2008-12-15 07:00:00.000               31     1.8
## 4 33042889055    true 2008-12-15 08:00:00.000               33     3.4
## 5 33042889056    true 2008-12-15 09:00:00.000               35     2.2
## 6 33042889057    true 2008-12-15 10:00:00.000               36     3.6
##   height.raw sensor.type individual.taxon.canonical.name tag.local.identifier
## 1     987.00         gps                    Panthera leo             GSM05751
## 2     983.90         gps                    Panthera leo             GSM05751
## 3     978.80         gps                    Panthera leo             GSM05751
## 4     990.30         gps                    Panthera leo             GSM05751
## 5     987.63         gps                    Panthera leo             GSM05751
## 6     987.52         gps                    Panthera leo             GSM05751
##   individual.local.identifier
## 1                        1003
## 2                        1003
## 3                        1003
## 4                        1003
## 5                        1003
## 6                        1003
##                                                   study.name
## 1 African Lion - Panthera Lion Central Kalahari Game Reserve
## 2 African Lion - Panthera Lion Central Kalahari Game Reserve
## 3 African Lion - Panthera Lion Central Kalahari Game Reserve
## 4 African Lion - Panthera Lion Central Kalahari Game Reserve
## 5 African Lion - Panthera Lion Central Kalahari Game Reserve
## 6 African Lion - Panthera Lion Central Kalahari Game Reserve
##                tstamp       date year month day hour min sec   idyear
## 1 2008-12-15 05:00:00 2008-12-15 2008    12  15    5   0   0 10032008
## 2 2008-12-15 06:00:00 2008-12-15 2008    12  15    6   0   0 10032008
## 3 2008-12-15 07:00:00 2008-12-15 2008    12  15    7   0   0 10032008
## 4 2008-12-15 08:00:00 2008-12-15 2008    12  15    8   0   0 10032008
## 5 2008-12-15 09:00:00 2008-12-15 2008    12  15    9   0   0 10032008
## 6 2008-12-15 10:00:00 2008-12-15 2008    12  15   10   0   0 10032008
##                 POSIX       X        Y
## 1 2008-12-15 05:00:00 2891248 -7700646
## 2 2008-12-15 06:00:00 2892276 -7699642
## 3 2008-12-15 07:00:00 2892975 -7699222
## 4 2008-12-15 08:00:00 2892986 -7699227
## 5 2008-12-15 09:00:00 2893087 -7699449
## 6 2008-12-15 10:00:00 2893092 -7699448

# First turn the SF into a data frame
animal_cleanedDF <- as.data.frame(animalSF_proj)

# Now save this into a csv file, giving the name of the projection in the filename
write.csv(animal_cleanedDF, 'animal_CleanData_projected_ESRI102014.csv')

# Create a set of trajectories, one per each individual
animal_ltraj <- as.ltraj(xy=animalSF_proj[,c('X','Y')],        # spatial Coordinates
                       date = animalSF_proj$POSIX,           # timestamp in POSIX format, including date/time
                       
                       id = animalSF_proj$individual.local.identifier)                # individual IDs

# Let's plot these trajectories to see if it worked
plot(animal_ltraj)

# Let's focus on 10 African Lions for this analysis, I will choose the individual 1001 to 1010. Include the the long-distance traveller and sedentary.
animal1 <- animal_ltraj[[1]]
animal2 <- animal_ltraj[[2]]
animal3 <- animal_ltraj[[3]]
animal4 <- animal_ltraj[[4]]
animal5 <- animal_ltraj[[5]]
animal6 <- animal_ltraj[[6]]
animal7 <- animal_ltraj[[7]]
animal8 <- animal_ltraj[[8]]
animal9 <- animal_ltraj[[9]]
animal10 <- animal_ltraj[[10]]
animal_list <- list(animal1, animal2, animal3, animal4, animal5,
                    animal6, animal7, animal8, animal9, animal10)

Figure of the duration, length and speed (choosen 10)

for (i in 1:length(animal_list)) {
  animal <- animal_list[[i]]
  animal$speedMS <- animal$dist / animal$dt
  animal$speedKmH <- animal$speedMS * 3.6
  
  duration_plot <- ggplot(animal, aes(x=dt)) +
    geom_histogram() +
    theme_bw() +
    xlab("Segment duration (s)") +
    ylab("Count") +
    ggtitle(paste("Segment Duration - Animal", i))
  
  length_plot <- ggplot(animal, aes(x=dist)) +
    geom_histogram() +
    theme_bw() +
    xlab("Segment length (m)") +
    ylab("Count") +
    ggtitle(paste("Segment Length - Animal", i))
  
  speed_plot <- ggplot(animal, aes(x=speedKmH)) +
    geom_histogram() +
    theme_bw() +
    xlab("Speed (km/h)") +
    ylab("Count") +
    ggtitle(paste("Speed - Animal", i))
  
  grid.newpage()
  pushViewport(viewport(layout = grid.layout(3, 1)))
  print(duration_plot, vp = viewport(layout.pos.row = 1, layout.pos.col = 1))
  print(length_plot, vp = viewport(layout.pos.row = 2, layout.pos.col = 1))
  print(speed_plot, vp = viewport(layout.pos.row = 3, layout.pos.col = 1))
}

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 2 rows containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 2 rows containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

cleaned_list <- list()

for (i in 1:length(animal_list)) {
  
  animal <- animal_list[[i]]
  
  before_plot <- ggplot(animal, aes(x=dt)) + 
    geom_boxplot() + 
    ggtitle(paste0("Before cleaning: animal", i)) + 
    theme_minimal()
  
  animal_cleaned <- animal[-which(animal$dt > 10000), ]
  cleaned_list[[i]] <- animal_cleaned
  
  after_plot <- ggplot(animal_cleaned, aes(x=dt)) + 
    geom_boxplot() + 
    ggtitle(paste0("After cleaning: animal", i)) + 
    theme_minimal()
  grid.arrange(before_plot, after_plot, nrow = 2)
}

## Warning: Removed 2 rows containing non-finite outside the scale range
## (`stat_boxplot()`).
## Removed 2 rows containing non-finite outside the scale range
## (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).
## Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).
## Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).
## Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).
## Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).
## Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).
## Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).
## Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).
## Removed 1 row containing non-finite outside the scale range (`stat_boxplot()`).

n <- length(animal_list)  # Number of individuals

# Initialize an empty data frame with columns for the statistics we want to calculate
results <- data.frame(
  dt_mean = numeric(n),
  dt_std = numeric(n),
  dist_mean = numeric(n),
  dist_std = numeric(n),
  speedKmH_mean = numeric(n),
  speedKmH_std = numeric(n)
)

# Step 3: For loop to calculate all the statistics

for (i in 1:length(animal_list)) {

  # Here I have copied the calculation code from step 1, but edited to include the index i
  
  # First we need to assign the data to a generic variable we will use in the for loop
  thisanimal <- animal_ltraj[[i]] # note how this is the same as above, except I replaced 1 with the index i

  # We also need to calculate speed, because we haven't done that for all ducks yet
  thisanimal$speedKmH <- (thisanimal$dist/thisanimal$dt)*3.6

  # Then we need to get rid of all NA values in the respective columns, using the two concepts you already know, is.na and | operator
  dtNAs <- which(is.na(thisanimal$dt)==TRUE | is.na(thisanimal$dist)==TRUE | is.na(thisanimal$speedKmH)==TRUE)
  # We remove rows with NAs
  thisanimal <- thisanimal[-dtNAs,]

  # Finally we calculate statistics for this duck, but unlike above, we now put these into the results data frame, in row i
  results$dt_mean[i] <- mean(thisanimal$dt)
  results$dt_std[i] <- sd(thisanimal$dt)
  results$dist_mean[i] <- mean(thisanimal$dist)
  results$dist_std[i] <- sd(thisanimal$dist)
  results$speedKmH_mean[i] <- mean(thisanimal$speedKmH)
  results$speedKmH_std[i] <- sd(thisanimal$speedKmH)

} 

# We don't have to initialise the results data frame here, because we already have it. But we will add two columns to it, like this:
results$R2n_mean <- NA
results$R2n_std <- NA

# Now let's use a for loop to populate these two columns. The loop is very similar to the previous one:

for (i in 1:length(animal_list)) {

  # Assign data to the generic variable
  thisanimal <- animal_ltraj[[i]] 

  # Then we need to get rid of all NA values in the respective column
  R2nNAs <- which(is.na(thisanimal$R2n)==TRUE)
  
  # Here comes the new If statement: if there are some NAs, we remove them, if not, we don't do anything.
  if (length(R2nNAs)>0) {
    thisanimal <- thisanimal[-R2nNAs,]
  }

  # Finally we calculate statistics for this duck, but unlike above, we now put these into the results data frame, in row i
  results$R2n_mean[i] <- mean(thisanimal$R2n)
  results$R2n_std[i] <- sd(thisanimal$R2n)

} 

4. Net Squared Displacement (NSD)

for (i in 1:length(animal_list)) {
  # Create the plot for each animal
  plot <- ggplot(animal_ltraj[[i]], aes(x=date, y=R2n)) + 
    geom_point() + 
    theme_bw() + 
    xlab("Date") + 
    ylab("R2n") + 
    ggtitle(paste("Net square displacement over time for animal", i))
  
  # Print the plot for each animal
  print(plot)
}

## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_point()`).

5. Turning Angle and Heading

# Take one duck
thisanimal <- animal_ltraj[[1]]

# Convert relative angle from radians into degrees
thisanimal$turningAngle <- thisanimal$rel.angle * 180 / pi
# Check it
head(thisanimal)

##            x        y                date          dx          dy       dist
## 1990 2703787 -7749541 2009-12-07 19:30:00  29.3199238  -56.009937  63.220020
## 1991 2703816 -7749597 2009-12-07 20:00:00  -6.4679256   -1.390088   6.615618
## 1992 2703810 -7749598 2009-12-07 20:30:00   0.3210739   -1.791937   1.820475
## 1993 2703810 -7749600 2009-12-07 21:00:00  -0.1821387    3.435114   3.439939
## 1994 2703810 -7749597 2009-12-07 21:30:00 -27.5632814   44.996767  52.767827
## 1995 2703782 -7749552 2009-12-07 22:00:00 112.9555376 -230.134102 256.360407
##        dt       R2n abs.angle  rel.angle turningAngle
## 1990 1800    0.0000 -1.088544         NA           NA
## 1991 1800 3996.7709 -2.929893 -1.8413489   -105.50152
## 1992 1800 3816.9766 -1.393501  1.5363920     88.02878
## 1993 1800 4040.6796  1.623769  3.0172700    172.87684
## 1994 1800 3637.4091  2.120401  0.4966316     28.45490
## 1995 1800  136.6857 -1.114516  3.0482686    174.65293

# Build a circular plot
angleplot <- ggplot(thisanimal, aes(x=turningAngle)) + geom_histogram(binwidth=2)+coord_polar(start = pi)+
  theme_bw()+scale_x_continuous(limits = c(-180,180),breaks=c(-180,-90,0,90,180))+xlab("Turning Angle")+ylab("Count")+
  ggtitle("Turning angle for animal 1")

# Export as a figure
png(filename="animal1_TurningAngle.png")
print(angleplot)

## Warning: Removed 19 rows containing non-finite outside the scale range
## (`stat_bin()`).

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

dev.off()

## quartz_off_screen 
##                 2

# Remove NAs 
turnNAs <- which(is.na(thisanimal$turningAngle)==TRUE)
  
# Here comes the new If statement: if there are some NAs, we remove them, if not, we don't do anything.
if (length(turnNAs)>0) {
    thisanimal <- thisanimal[-turnNAs,]
}

# Calculate statistics - because angles are circular, we need to use circular mean and circular std, from package circular
turningAngle_mean <- mean(circular(thisanimal$turningAngle))
turningAngle_std <- sd(circular(thisanimal$turningAngle))

# Take a lion
thisanimal <- animal_ltraj[[2]]

# Check to see in which columns we have the x and y coordinates
head(thisanimal)

##             x        y                date          dx         dy       dist
## 80158 2870833 -7710273 2009-08-14 18:00:00   1.9234920  1.2928059  2.3175781
## 80159 2870835 -7710271 2009-08-14 19:00:00  -0.0852056  0.4467048  0.4547584
## 80160 2870835 -7710271 2009-08-14 20:00:00 -14.5114015 -9.2493346 17.2084562
## 80161 2870820 -7710280 2009-08-14 21:00:00   8.9965413  4.0308027  9.8582516
## 80162 2870829 -7710276 2009-08-14 22:00:00  -1.5971638  1.5471826  2.2236695
## 80163 2870828 -7710275 2009-08-14 23:00:00   2.4492185 -6.0142335  6.4938183
##         dt        R2n  abs.angle rel.angle
## 80158 3600   0.000000  0.5917644        NA
## 80159 3600   5.371168  1.7592748  1.167510
## 80160 3600   6.405194 -2.5741375  1.949773
## 80161 3600 217.005304  0.4212221  2.995360
## 80162 3600  25.620784  2.3720888  1.950867
## 80163 3600  31.544310 -1.1840668  2.727030

# Ok, they are called x and y. 

# So now we will take a loop that will go through this trajecotry. We need to know the number of points in the trajectory,
# we can just check the length of one of the coordinates.
nPoints <- length(thisanimal$x)

# Also we need to add a new field for the heading to thisduck data frame
thisanimal$heading <- NA

# For loop through the trajectory, but not to the end, we skip the final point!

for (i in 1:(nPoints-1)) {
  
  # Let's read the coordinates from two consecutive points, i and i+1
  x1 <- thisanimal$x[i]
  x2 <- thisanimal$x[i+1]
  y1 <- thisanimal$y[i]
  y2 <- thisanimal$y[i+1]
  
  # And calculate the heading as arctan of the two values, plus convert to degrees
  thisanimal$heading[i] <- atan2((x2-x1),(y2-y1)) * 180 / pi
  
}

# Build a circular plot
angleplot <- ggplot(thisanimal, aes(x=heading)) + geom_histogram(binwidth=2)+coord_polar(start = pi)+
  theme_bw()+scale_x_continuous(limits = c(-180,180),breaks=c(-180,-90,0,90,180))+
  xlab("Heading")+ylab("Count")+ggtitle("Heading for duck 2")

# Export as a figure
png(filename="Animal2_Heading.png")
print(angleplot)

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

dev.off()

## quartz_off_screen 
##                 2

# Remove NAs - there will be at least one, for the last point through which we didn't loop
headingNAs <- which(is.na(thisanimal$heading)==TRUE)
  
# If there are some NAs, we remove them, if not, we don't do anything.
if (length(headingNAs)>0) {
    thisanimal <- thisanimal[-headingNAs,]
}

# Calculate statistics - because angles are circular, we need to use circular mean and circular std, from package circular
heading_mean <- mean(circular(thisanimal$heading))
heading_std <- sd(circular(thisanimal$heading))

turning_summary <- data.frame(
  animal_id = numeric(),
  mean_angle = numeric(),
  std_angle = numeric()
)

for (i in 1:length(animal_list)) {
  animal <- animal_list[[i]]
  
  if (!all(c("x", "y") %in% names(animal))) next
  animal$turning <- NA
  
  for (j in 2:(nrow(animal) - 1)) {
    x1 <- animal$x[j - 1]
    x2 <- animal$x[j]
    x3 <- animal$x[j + 1]
    y1 <- animal$y[j - 1]
    y2 <- animal$y[j]
    y3 <- animal$y[j + 1]
    
    angle1 <- atan2(y2 - y1, x2 - x1)
    angle2 <- atan2(y3 - y2, x3 - x2)
    turn <- angle2 - angle1
    
    if (turn > pi) turn <- turn - 2 * pi
    if (turn < -pi) turn <- turn + 2 * pi
    
    animal$turning[j] <- turn * 180 / pi
  }
  
  animal_cleaned <- animal[!is.na(animal$turning), ]
  
  mean_turn <- mean(circular(animal_cleaned$turning, units = "degrees"))
  sd_turn <- sd(circular(animal_cleaned$turning, units = "degrees"))
  
  turning_summary <- rbind(turning_summary, data.frame(
    animal_id = i,
    mean_angle = mean_turn,
    std_angle = sd_turn
  ))

  angleplot <- ggplot(animal_cleaned, aes(x = turning)) +
    geom_histogram(binwidth = 2) +
    coord_polar(start = pi) +
    theme_bw() +
    scale_x_continuous(limits = c(-180, 180), breaks = c(-180, -90, 0, 90, 180)) +
    xlab("Turning Angle") + ylab("Count") +
    ggtitle(paste("Turning Angle for Animal", i))
  
  print(angleplot)
  cat(paste0("Animal ", i, " – Mean Turning Angle: ", round(mean_turn, 2),
             "°, Std Dev: ", round(sd_turn, 2), "°\n\n"))
}

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 1 – Mean Turning Angle: -0.29°, Std Dev: 1.93°

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 2 – Mean Turning Angle: -0.6°, Std Dev: 2.03°

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 3 – Mean Turning Angle: 2.99°, Std Dev: 1.53°

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 4 – Mean Turning Angle: 25.42°, Std Dev: 2.45°

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 5 – Mean Turning Angle: -131.31°, Std Dev: 3.13°

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 6 – Mean Turning Angle: 37.55°, Std Dev: 1.51°

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 7 – Mean Turning Angle: 20.84°, Std Dev: 2.62°

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 8 – Mean Turning Angle: 123.15°, Std Dev: 2.71°

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 9 – Mean Turning Angle: -175.56°, Std Dev: 1.94°

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Animal 10 – Mean Turning Angle: 175.32°, Std Dev: 2.34°

# We don't have to initialise the results data frame here, because we already have it. But we will add two columns to it, like this:
results$heading_mean <- NA
results$heading_std <- NA

# Here start the nested loops

# External For loop for all the ducks - note the index is now called j!
for (j in 1:length(animal_list)) {

  # Take a lion
  thisanimal <- animal_ltraj[[j]]

  # So now we will take a loop that will go through this trajecotry. We need to know the number of points in the   trajectory,
  # we can just check the length of one of the coordinates.
  nPoints <- length(thisanimal$x)

  # Also we need to add a new field for the heading to thisanimal data frame
  thisanimal$heading <- NA

  # Internal for loop through the trajectory, but not to the end, we skip the final point!
  # This loop has index i
  for (i in 1:(nPoints-1)) {
  
    # Let's read the coordinates from two consecutive points, i and i+1
    x1 <- thisanimal$x[i]
    x2 <- thisanimal$x[i+1]
    y1 <- thisanimal$y[i]
    y2 <- thisanimal$y[i+1]
  
    # And calculate the heading as arctan of the two values, plus convert to degrees
    thisanimal$heading[i] <- atan2((x2-x1),(y2-y1)) * 180 / pi
  
  } # end for i

  # Build a circular plot
  angleplot <- ggplot(thisanimal, aes(x=heading)) + geom_histogram(binwidth=2)+coord_polar(start = pi)+
    theme_bw()+scale_x_continuous(limits = c(-180,180),breaks=c(-180,-90,0,90,180))+
    xlab("Heading")+ylab("Count")+ggtitle(paste("Heading for animal ",j,sep=""))
  
  # Export as a figure 
  #png(filename=paste("Animal",j,"_Heading.png",sep=""))
  print(angleplot)
  #dev.off()

  # Remove NAs - there will be at least one, for the last point through which we didn't loop
  headingNAs <- which(is.na(thisanimal$heading)==TRUE)
  
  # Here comes the new If statement: if there are some NAs, we remove them, if not, we don't do anything.
  if (length(headingNAs)>0) {
    thisanimal <- thisanimal[-headingNAs,]
  }

  # Calculate statistics - because angles are circular, we need to use circular mean and circular std
  results$heading_mean[j] <- mean(circular(thisanimal$heading))
  results$heading_std[j] <- sd(circular(thisanimal$heading))
  
} # end for j

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_bin()`).

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_bin()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_bin()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_bin()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_bin()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_bin()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_bin()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_bin()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_bin()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_bin()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_bar()`).

# Check results
results

##       dt_mean      dt_std  dist_mean  dist_std speedKmH_mean speedKmH_std
## 1    2864.155    4921.787   736.2872  1270.044     1.3961524    2.3914207
## 2    7258.505   34250.186   962.3112  2922.488     0.9232217    1.6871668
## 3    6477.465   12728.197  1404.3166  2178.671     1.2778381    1.8743605
## 4   98387.025  330931.434  7274.8486 10273.839     0.3779828    0.6003488
## 5   10479.860   16850.670  1561.3494  4530.184     0.4102659    0.9691038
## 6  481693.043 1159510.388 12477.7770 22807.435     0.5326418    0.8773748
## 7   66831.111  174078.629  5539.7002  8103.403     0.4435427    0.7994241
## 8   28078.062  103455.111  2611.2533  6216.715     0.3440738    0.7046752
## 9   31443.061   96593.286  1408.3242  3627.827     0.2664249    0.7473899
## 10  19825.389   52526.566  2543.1237  6588.879     0.5157002    1.0516710
##      R2n_mean    R2n_std heading_mean heading_std
## 1  1246054400 2182639514   2.11934272    2.875035
## 2  1532098517 2013510350  -1.76883741    3.148480
## 3  1310816225 1097295839  -0.35291785    2.997561
## 4  3690232022 2715599705   0.90269072    2.693817
## 5  1082938101 1093634354   1.13997808    2.661001
## 6  1128149017 1169569646  -1.77421529    1.806196
## 7   450897353  546744742   2.27406183    2.813500
## 8   902467841  943075197  -0.01773315    2.439466
## 9   549330014  284666991   3.13593078    2.539468
## 10 1891362359 2542130782  -0.34558583    2.180880

6. Radius of Gyration

# Prepare the field in the results data frame
results$gyration <- NA

# For loop across all ducks
for (i in 1:length(animal_list)) {

  # Pick a duck
  thisanimal <- animal_ltraj[[i]]

  # Find number of points in this trajectory
  noPoints <- length(thisanimal$x)

  # Calculate the mean centre of all locations
  centreX <- sum(thisanimal$x)/noPoints
  centreY <- sum(thisanimal$y)/noPoints

  # Distance of each point to (centreX, centreY)
  thisanimal$distTocentre <- sqrt((centreX-thisanimal$x)^2+(centreY-thisanimal$y)^2)

  # Calculate radius of gyration and insert it into the results data frame
  results$gyration[i] <- sqrt(sum(thisanimal$distTocentre)/noPoints)
  
}
# Check results
results

##       dt_mean      dt_std  dist_mean  dist_std speedKmH_mean speedKmH_std
## 1    2864.155    4921.787   736.2872  1270.044     1.3961524    2.3914207
## 2    7258.505   34250.186   962.3112  2922.488     0.9232217    1.6871668
## 3    6477.465   12728.197  1404.3166  2178.671     1.2778381    1.8743605
## 4   98387.025  330931.434  7274.8486 10273.839     0.3779828    0.6003488
## 5   10479.860   16850.670  1561.3494  4530.184     0.4102659    0.9691038
## 6  481693.043 1159510.388 12477.7770 22807.435     0.5326418    0.8773748
## 7   66831.111  174078.629  5539.7002  8103.403     0.4435427    0.7994241
## 8   28078.062  103455.111  2611.2533  6216.715     0.3440738    0.7046752
## 9   31443.061   96593.286  1408.3242  3627.827     0.2664249    0.7473899
## 10  19825.389   52526.566  2543.1237  6588.879     0.5157002    1.0516710
##      R2n_mean    R2n_std heading_mean heading_std gyration
## 1  1246054400 2182639514   2.11934272    2.875035 147.7880
## 2  1532098517 2013510350  -1.76883741    3.148480 163.5427
## 3  1310816225 1097295839  -0.35291785    2.997561 139.6413
## 4  3690232022 2715599705   0.90269072    2.693817 188.2696
## 5  1082938101 1093634354   1.13997808    2.661001 145.4142
## 6  1128149017 1169569646  -1.77421529    1.806196 155.8137
## 7   450897353  546744742   2.27406183    2.813500 131.1873
## 8   902467841  943075197  -0.01773315    2.439466 159.2748
## 9   549330014  284666991   3.13593078    2.539468  95.2824
## 10 1891362359 2542130782  -0.34558583    2.180880 184.4834

# This is the last thing we will add to results, so let's also export these into a file.
write.csv(results, "Animal_results.csv", row.names=FALSE)

Analytical Interpretation

What time period do we have the data from?

The dataset spans from 2008 to 2012, covering a multi-year period across multiple seasons. This range provides a robust foundation to explore both short-term behavioural rhythms and long-term patterns. Each individual has a slightly different data coverage, but the full temporal extent allows for comparative and time-based movement analytics.

What can you say about individuals’ movement patterns based on these results?

Across the 10 individuals, movement behaviour diverges significantly in both spatial extent and temporal regularity. Based on net squared displacement (NSD), radius of gyration, trajectory maps, and heading distributions, we identify three broad behavioural categories:

1. long distant travel individulas (e.g., Animals 1, 2, 4, 10): These individuals demonstrate large NSD values, high radius of gyration, and sprawling trajectories. Their movements extend over tens of kilometers, suggesting exploratory or migratory behaviour, likely in search of resources or seasonal habitats.
2. Sedentary individuals (e.g., Animals 5, 6, 9): These exhibit minimal displacement from origin, low gyration, and highly clustered movement paths. Their segment durations are variable, but their speeds are generally low (below 6 km/h). These patterns point to site fidelity, residency, or constrained movement within a home range.
3. Cyclic individulas (notably Animal 3): This individual displays a recurrent NSD waveform—clear cycles of movement away from and return to the origin point. This may reflect commuting behaviour, patrol activity, or seasonal site returns, and is supported by moderately concentrated heading values.

Are there any individuals that are behaving differently than the others?

Yes—several individuals stand out for their unique or potentially anomalous movement signatures:

1. Animal 2 has some of the highest values for displacement and speed, with a mean segment duration over 7000 seconds and maximum values reaching close to 1 million seconds. The corresponding heading plot is nearly uniform, which could imply non-directed dispersal. Alternatively, this might reflect tracking inconsistencies or data irregularities, as the duration and spread are far beyond most others.
2. Animal 6 has sparse data points, but high mean segment duration and large individual displacements, yet very low speeds. This combination may indicate data gaps, GPS error, or inactivity. Although it shows large distances moved, the low speeds and poor resolution reduce interpretability.
3. Animal 3, as previously mentioned, exhibits distinct periodicity in its NSD and movement vectors. This behaviour contrasts with more monotonic dispersal or static patterns seen in others.

How do turning angle and heading correspond to the map of individuals’ trajectories?

1. Animals 1, 2, and 3 have turning angle distributions tightly centered around 0°, with low standard deviations (less than 2°). This indicates mostly straight-line movement with minimal turning, which aligns with their relatively linear trajectories on the map. Their heading plots further support this by showing concentrated directions, especially for Animal 3.
2. Animals 4 and 7 have moderate mean turning angles (20°–25°) with higher spread, suggesting frequent directional adjustments, likely representing patrol-like or area-restricted search behaviour within defined territories.
3. Animals 5, 8, 9, and 10 exhibit extreme turning angles (e.g., -175°, +175°, ±130°), meaning they often reverse or loop back sharply. This behaviour corresponds to highly compact or criss-crossed trajectories seen in the spatial plots. For example, Animal 5’s turning angle of -131° suggests a strong bias toward leftward circular movement, which may be due to environmental boundaries.
4. Animal 6, while having a large mean turning angle (≈37°), suffers from limited data and low angular resolution. This outlier suggests either sparse tracking or highly constrained movement that doesn’t allow diverse angles to emerge.

Together, the turning angle plots reflect the tight vs. straight curves, while the heading plots indicate overall movement direction. For example, Animals 9 and 10 both show turning means close to ±180°, indicating repeated sharp reversals, which matches their dense, non-linear paths on the trajectory map.

1. Tightly grouped angles and consistent headings reflect purposeful, often long-range movement;
2. Broad or bimodal turning patterns correspond to local searching or looping;
3. Extremely sharp angles suggest environmental barriers, territorial boundaries, or artifacts of limited space.

In short, integrating turning angle and heading statistics with trajectory maps greatly enhances our ability to infer movement modes—whether directed travel, exploratory behaviour, or site fidelity.

Compare results of net squared displacement and radius of gyration with maps of all individuals

These two spatial dispersion metrics offer distinct yet complementary insights:

1. High NSD and high gyration (e.g., Animals 1, 2, 4, 10): These animals traveled far from their starting locations and exhibited wide-ranging activity around new centers. The trajectory maps support this, showing dispersed locations across large areas.
2. High NSD, low gyration (e.g., Animal 3): This indicates linear or directional movement away from origin, with minimal lateral spread. The cyclical pattern in NSD supports this interpretation—straight but returning movements.
3. Low NSD and gyration (e.g., Animals 5, 6, 9): These animals remained close to their origin and had tightly bounded paths. In maps, their trajectories appear as dense point clusters, consistent with localized foraging or denning behaviour.

By aligning spatial metrics with visual maps, we can validate movement interpretations and differentiate functional movement strategies—e.g., exploration, migration, residency.

Additional Reflections

This analysis assumes that the GPS data accurately reflects real movement, but data irregularities—such as long time gaps (Animal 6) or abnormally long segments (Animal 2)—can bias estimates of speed and displacement. Cleaning these outliers was essential, but limitations remain.

Future work could combine this movement data with environmental layers to model habitat preferences or movement constraints in a landscape ecology context.

Reference

Bastille-Rousseau, G., Potts, J.R., Yackulic, C.B. et al. Flexible characterization of animal movement pattern using net squared displacement and a latent state model. Mov Ecol 4, 15 (2016). https://doi.org/10.1186/s40462-016-0080-y

Codling, E.A., Plank, M.J., & Benhamou, S. (2008). Random walk models in biology. Journal of the Royal Society Interface, 5(25), 813–834. https://doi.org/10.1098/rsif.2008.0014

Nathan, R. et al. (2008). A movement ecology paradigm for unifying organismal movement research. Proceedings of the National Academy of Sciences, 105(49), 19052–19059. https://doi.org/10.1073/pnas.0800375105